# Cleaning the Work Environment ####
rm(list=ls(all=T))
gc()

getwd()

library(sjPlot)
library(multcomp)
library(stargazer)

#### Loading the Dataset ####

s1 <- readRDS(file = "SEEBSS - Study 1.rds")

#### Power Analysis #####

pwr::pwr.t.test(d=0.30, sig.level=.05, power = .80, type = 'two.sample')

#### Data Preperation ####

table(s1$Q31_1, exclude = NULL)
s1$Q31_1[s1$Q31_1==99] <- NA
names(s1)[names(s1)=="Q31_1"] <- "Self_Conservatism"
table(s1$Self_Conservatism, exclude = NULL)

table(s1$Q36, exclude = NULL)
s1$Q36[s1$Q36==99] <- NA
s1$Q36[s1$Q36==2] <- 0
names(s1)[names(s1)=="Q36"] <- "sex" # 1 - female, 0 - male
table(s1$sex, exclude = NULL)
s1$sex <- 1 + s1$sex*-1 # recoding # 0 - female, 1 - male
table(s1$sex, exclude = NULL) 

s1$Q38[s1$Q38==99] <- NA
names(s1)[names(s1)=="Q38"] <- "education"
table(s1$education, exclude = NULL)

s1$Q7_1[s1$Q7_1==99] <- NA
names(s1)[names(s1)=="Q7_1"] <- "self_left_right"
table(s1$self_left_right, exclude = NULL)

s1$Affiliation <- NA
s1$Affiliation[s1$Q3==99] <- 0
s1$Affiliation[s1$Q3==1] <- 2
s1$Affiliation[s1$Q3==2] <- 1
table(s1$Affiliation, exclude = NULL)

table(s1$Q4_8_TEXT, exclude = NULL)

table(s1$Q4, exclude = NULL)
s1$Partisanship <- NA
s1$Partisanship[s1$Q4==1] <- "Justice and Development Party (AKP)"
s1$Partisanship[s1$Q4==2] <- "Republican People's Party (CHP)"
s1$Partisanship[s1$Q4==3] <- "Other Parties on the Left"
s1$Partisanship[s1$Q4==4] <- "Other Parties on the Right"
s1$Partisanship[s1$Q4==5] <- "Other Parties on the Right"
s1$Partisanship[s1$Q4==6] <- "Other Parties on the Right"
s1$Partisanship[s1$Q4==7] <- "Other Parties on the Left"
s1$Partisanship[s1$Q4_8_TEXT == "EMEP" | 
                  s1$Q4_8_TEXT == "Özgürlük ve dayanışma partisi" | 
                  s1$Q4_8_TEXT == "TKP" | 
                  s1$Q4_8_TEXT == "Türkiye işçi partisi" | 
                  s1$Q4_8_TEXT == "Türkiye Komünist Partisi (TKP)" |
                  s1$Q4_8_TEXT == "Türkiye Komün"] <- "Other Parties on the Left"
s1$Partisanship[s1$Q4_8_TEXT == "Ldp" | s1$Q4_8_TEXT == "Hepar" |
                  s1$Q4_8_TEXT == "Hüdapar" |
                  s1$Q4_8_TEXT == "Ldp"] <- "Other Parties on the Right"
s1$Partisanship <- ifelse(is.na(s1$Partisanship),"Bipartisan or Non-Responsive", s1$Partisanship)
table(s1$Partisanship, exclude=NULL)
s1$Partisanship <- factor(s1$Partisanship, ordered = F) 
table(s1$Partisanship, exclude=NULL)

s1$bipartisan_na <- NA
s1$bipartisan_na[s1$Partisanship=="Bipartisan or Non-Responsive"] <- 1
s1$bipartisan_na[s1$Partisanship!="Bipartisan or Non-Responsive"] <- 0
table(s1$bipartisan_na, exclude = NULL)

s1$akp <- NA
s1$akp[s1$Partisanship=="Justice and Development Party (AKP)"] <- 1
s1$akp[s1$Partisanship!="Justice and Development Party (AKP)"] <- 0
table(s1$akp, exclude = NULL)

s1$chp <- NA
s1$chp[s1$Partisanship=="Republican People's Party (CHP)"] <- 1
s1$chp[s1$Partisanship!="Republican People's Party (CHP)"] <- 0
table(s1$chp, exclude = NULL)

s1$other_right <- NA
s1$other_right[s1$Partisanship=="Other Parties on the Right"] <- 1
s1$other_right[s1$Partisanship!="Other Parties on the Right"] <- 0
table(s1$other_right, exclude = NULL)

s1$other_left <- NA
s1$other_left[s1$Partisanship=="Other Parties on the Left"] <- 1
s1$other_left[s1$Partisanship!="Other Parties on the Left"] <- 0
table(s1$other_left, exclude = NULL)

mean(s1$akp, na.rm = T)
s1$beta_akp <- (s1$akp - mean(s1$akp, na.rm = T))
table(s1$beta_akp, exclude = NULL)

mean(s1$chp, na.rm = T)
s1$beta_chp <- (s1$chp - mean(s1$chp, na.rm = T))
table(s1$beta_chp, exclude = NULL)

mean(s1$other_right, na.rm = T)
s1$beta_other_right <- (s1$other_right - mean(s1$other_right, na.rm = T))
table(s1$beta_other_right, exclude = NULL)

mean(s1$other_left, na.rm = T)
s1$beta_other_left <- (s1$other_left - mean(s1$other_left, na.rm = T))
table(s1$beta_other_left, exclude = NULL)

table(s1$Group, exclude = NULL)

s1$Cancer_Prime <- NA
s1$Cancer_Prime[s1$Group == "Treatment 1"] <- 1
s1$Cancer_Prime <- ifelse(is.na(s1$Cancer_Prime), 0, s1$Cancer_Prime)
table(s1$Cancer_Prime, exclude = NULL)

s1$Terror_Prime <- NA
s1$Terror_Prime[s1$Group == "Treatment 2"] <- 1
s1$Terror_Prime <- ifelse(is.na(s1$Terror_Prime), 0, s1$Terror_Prime)
table(s1$Terror_Prime, exclude = NULL)

s1$Domestic_Violence_Prime <- NA
s1$Domestic_Violence_Prime[s1$Group == "Treatment 3"] <- 1
s1$Domestic_Violence_Prime <- ifelse(is.na(s1$Domestic_Violence_Prime), 0, s1$Domestic_Violence_Prime)
table(s1$Domestic_Violence_Prime, exclude = NULL)

table(s1$Q14_4, exclude = NULL)
s1$Q14_4 <- ifelse(s1$Q14_4>11, NA, s1$Q14_4)
s1$Q14_4[s1$Group=="Control"] <- 0

table(s1$Q14_5, exclude = NULL)
s1$Q14_5 <- ifelse(s1$Q14_5>11, NA, s1$Q14_5)
s1$Q14_5[s1$Group=="Control"] <- 0

table(s1$Q14_6, exclude = NULL)
s1$Q14_6 <- ifelse(s1$Q14_6>11, NA, s1$Q14_6)
s1$Q14_6[s1$Group=="Control"] <- 0

s1$cancer_prime_fear <- s1$Q14_4 + s1$Q14_5 + s1$Q14_6
table(s1$cancer_prime_fear, exclude = NULL)

table(s1$Q15_4, exclude = NULL)
s1$Q15_4 <- ifelse(s1$Q15_4>11, NA, s1$Q15_4)
s1$Q15_4[s1$Group=="Control"] <- 0

table(s1$Q15_5, exclude = NULL)
s1$Q15_5 <- ifelse(s1$Q15_5>11, NA, s1$Q15_5)
s1$Q15_5[s1$Group=="Control"] <- 0

table(s1$Q15_6, exclude = NULL)
s1$Q15_6 <- ifelse(s1$Q15_6>11, NA, s1$Q15_6)
s1$Q15_6[s1$Group=="Control"] <- 0

s1$terror_prime_fear <- s1$Q15_4 + s1$Q15_5 + s1$Q15_6
table(s1$terror_prime_fear, exclude = NULL)

table(s1$Q16_4, exclude = NULL)
s1$Q16_4 <- ifelse(s1$Q16_4>11, NA, s1$Q16_4)
s1$Q16_4[s1$Group=="Control"] <- 0

table(s1$Q16_5, exclude = NULL)
s1$Q16_5 <- ifelse(s1$Q16_5>11, NA, s1$Q16_5)
s1$Q16_5[s1$Group=="Control"] <- 0

table(s1$Q16_6, exclude = NULL)
s1$Q16_6 <- ifelse(s1$Q16_6>11, NA, s1$Q16_6)
s1$Q16_6[s1$Group=="Control"] <- 0

s1$domestic_violence_prime_fear <- s1$Q16_4 + s1$Q16_5 + s1$Q16_6
table(s1$domestic_violence_prime_fear, exclude = NULL)

s1$priming_fear <- rowSums(s1[, c("cancer_prime_fear", "terror_prime_fear", "domestic_violence_prime_fear")], na.rm = T)
table(s1$priming_fear, exclude = NULL)
s1$priming_fear <- ifelse(is.na(s1$cancer_prime_fear) & is.na(s1$terror_prime_fear) & is.na(s1$domestic_violence_prime_fear), NA, s1$priming_fear)

#### Manipulation Check ####

s1$Q17_1[s1$Q17_1==99] <- NA
s1$Q17_2[s1$Q17_2==99] <- NA

s1$MC1 <- NA
s1$MC1[s1$Q17_1==6 & s1$Q17_2==8] <- 1
s1$MC1[(s1$Q17_1!=6 | s1$Q17_2!=8) & (!is.na(s1$Q17_1) & !is.na(s1$Q17_2))] <- 0
s1$MC1 <- ifelse(is.na(s1$MC1), 0, s1$MC1)
table(s1$MC1, exclude = NULL)

s1$Q18_1[s1$Q18_1==99] <- NA
s1$Q18_2[s1$Q18_2==99] <- NA

s1$MC2 <- NA
s1$MC2[s1$Q18_1==1 & s1$Q18_2==4] <- 1
s1$MC2[(s1$Q18_1!=1 | s1$Q18_2!=4) & (!is.na(s1$Q18_1) & !is.na(s1$Q18_2))] <- 0
s1$MC2 <- ifelse(is.na(s1$MC2), 0, s1$MC2)
table(s1$MC2, exclude = NULL)

s1$Q19_1[s1$Q19_1==99] <- NA
s1$Q19_2[s1$Q19_2==99] <- NA

s1$MC3 <- NA
s1$MC3[s1$Q19_1==1 & s1$Q19_2==4] <- 1
s1$MC3[(s1$Q19_1!=1 | s1$Q19_2!=4) & (!is.na(s1$Q19_1) & !is.na(s1$Q19_2))] <- 0
s1$MC3 <- ifelse(is.na(s1$MC3), 0, s1$MC3)
table(s1$MC3, exclude = NULL)

s1$Manipulation_Check <- s1$MC1 + s1$MC2+ s1$MC3
table(s1$Manipulation_Check, s1$Group, exclude = NULL)

s1$Treatment1 <- NA
s1$Treatment1[s1$Group == "Treatment 1"] <- 1
s1$Treatment1[s1$Group == "Control"] <- 0
table(s1$Treatment1, exclude = NULL)

s1$Treatment2 <- NA
s1$Treatment2[s1$Group == "Treatment 2"] <- 1
s1$Treatment2[s1$Group == "Control"] <- 0
table(s1$Treatment2, exclude = NULL)

s1$Treatment3 <- NA
s1$Treatment3[s1$Group == "Treatment 3"] <- 1
s1$Treatment3[s1$Group == "Control"] <- 0
table(s1$Treatment3, exclude = NULL)

bartlett.test(s1$Manipulation_Check ~ s1$Treatment1)
t.test(s1$Manipulation_Check ~ s1$Treatment1, var.equal=FALSE, alternative="less")

bartlett.test(s1$Manipulation_Check ~ s1$Treatment2)
t.test(s1$Manipulation_Check ~ s1$Treatment2, var.equal=FALSE, alternative="less")

bartlett.test(s1$Manipulation_Check ~ s1$Treatment3)
t.test(s1$Manipulation_Check ~ s1$Treatment3, var.equal=TRUE, alternative="less")

#### Supplementary Appendix, Table 13. Study 1 descriptive statistics. ####

library(table1)

my.render.cont <- function(x) {
  with(stats.apply.rounding(stats.default(x), digits=2), c("",
                                                           "Mean (SD)"=sprintf("%s (&plusmn; %s)", MEAN, SD)))
}
my.render.cat <- function(x) {
  c("", sapply(stats.default(x), function(y) with(y,
                                                  sprintf("%d (%0.0f %%)", FREQ, PCT))))
}

table1(~ Self_Conservatism + priming_fear + sex + age + education + self_left_right | Group, data=s1, 
       render.continuous=my.render.cont, render.categorical=my.render.cat, overall=NULL)

treatment <- aov(Self_Conservatism ~ factor(Group), data = s1)
summary(treatment)

sex <- aov(sex ~ factor(Group), data = s1)
summary(sex)

age <- aov(age ~ factor(Group), data = s1)
summary(age)

education <- aov(education ~ factor(Group), data = s1)
summary(education)

self_left_right <- aov(self_left_right ~ factor(Group), data = s1)
summary(self_left_right)

#### Standardization ####

s1 <- subset(s1, subset = !is.na(Self_Conservatism))

mean(s1$Cancer_Prime, na.rm = T)
s1$beta_Cancer_Prime <- (s1$Cancer_Prime - mean(s1$Cancer_Prime, na.rm = T))
table(s1$beta_Cancer_Prime, exclude = NULL)

mean(s1$Terror_Prime, na.rm = T)
s1$beta_Terror_Prime <- (s1$Terror_Prime - mean(s1$Terror_Prime, na.rm = T))
table(s1$beta_Terror_Prime, exclude = NULL)

mean(s1$Domestic_Violence_Prime, na.rm = T)
s1$beta_Domestic_Violence_Prime <- (s1$Domestic_Violence_Prime - mean(s1$Domestic_Violence_Prime, na.rm = T))
table(s1$beta_Domestic_Violence_Prime, exclude = NULL)

mean(s1$sex, na.rm = T)
s1$beta_sex <- (s1$sex - mean(s1$sex, na.rm = T))
table(s1$beta_sex, exclude = NULL)

mean(s1$age, na.rm = T)
sd(s1$age, na.rm = T)
s1$beta_age <- (s1$age - mean(s1$age, na.rm = T))/(2*sd(s1$age, na.rm = T))
table(s1$beta_age, exclude = NULL)

mean(s1$education, na.rm = T)
sd(s1$education, na.rm = T)
s1$beta_education <- (s1$education - mean(s1$education, na.rm = T))/(2*sd(s1$education, na.rm = T))
table(s1$beta_education, exclude = NULL)

mean(s1$self_left_right, na.rm = T)
sd(s1$self_left_right, na.rm = T)
s1$beta_self_left_right <- (s1$self_left_right - mean(s1$self_left_right, na.rm = T))/(2*sd(s1$self_left_right, na.rm = T))
table(s1$beta_self_left_right, exclude = NULL)

mean(s1$Self_Conservatism, na.rm = T)
sd(s1$Self_Conservatism, na.rm = T)
s1$beta_Self_Conservatism <- (s1$Self_Conservatism - mean(s1$Self_Conservatism, na.rm = T))/(2*sd(s1$Self_Conservatism, na.rm = T))
table(s1$beta_Self_Conservatism, exclude = NULL)

mean(s1$priming_fear, na.rm = T)
sd(s1$priming_fear, na.rm = T)
s1$beta_priming_fear <- (s1$priming_fear - mean(s1$priming_fear, na.rm = T))/(2*sd(s1$priming_fear, na.rm = T))
table(s1$beta_priming_fear, exclude = NULL)

s1$region <- NA
s1$region[s1$State=='Istanbul'] <- "Istanbul"
s1$region[s1$State=='Tekirdag' |
            s1$State=='Edirne' |
            s1$State== 'Kirklareli' |
            s1$State== 'Balikesir' | 
            s1$State== 'Canakkale'] <- "West Marmara"
s1$region[s1$State=='Izmir' |
            s1$State== 'Aydin' | 
            s1$State==  'Denizli' |
            s1$State==  'Mugla' |
            s1$State==  'Manisa' |
            s1$State==  'Afyon' |
            s1$State== 'Kutahya' |
            s1$State== 'Usak'] <- "Aegean"
s1$region[s1$State=='Bursa' |
            s1$State== 'Eskisehir' |
            s1$State== 'Bilecik' |
            s1$State== 'Kocaeli' |
            s1$State== 'Sakarya' |
            s1$State==  'Duzce' |
            s1$State==  'Bolu' |
            s1$State== 'Yalova'] <- "East Marmara"
s1$region[s1$State=='Ankara' |
            s1$State==  'Konya' | s1$State== 'Karaman'] <- "West Anatolia"
s1$region[s1$State=='Antalya' |
            s1$State== 'Isparta' |
            s1$State==  'Adana' |
            s1$State==   'Mersin' |
            s1$State==   'Hatay' |
            s1$State==  'Kahramanmaras' |
            s1$State==   'Osmaniye' | s1$State== 'Burdur'] <- "Mediterranean"
s1$region[s1$State=='Kirikkale' |
            s1$State== 'Aksaray' |
            s1$State==  'Nigde' |
            s1$State==   'Nevsehir' |
            s1$State==   'Kirsehir' |
            s1$State==   'Kayseri' | 
            s1$State==  'Sivas' |
            s1$State==  'Yozgat'] <- "Central Anatolia"
s1$region[s1$State=='Zonguldak' |
            s1$State== 'Karabuk' |
            s1$State==  'Bartin' |
            s1$State==  'Kastamonu' |
            s1$State== 'Cankiri' |
            s1$State==  'Sinop' |
            s1$State==  'Samsun' |
            s1$State==  'Tokat' |
            s1$State==  'Corum' |
            s1$State== 'Amasya'] <- "West Black Sea"
s1$region[s1$State=='Ordu' |
            s1$State=='Trabzon' |
            s1$State=='Giresun' |
            s1$State== 'Rize' |
            s1$State== 'Artvin' |
            s1$State== 'Gumushane'] <- "East Black Sea"
s1$region[s1$State=='Erzurum' |
            s1$State=='Erzincan' |
            s1$State== 'Agri' |
            s1$State== 'Kars' | s1$State== 'Ardahan' | s1$State== 'Bayburt' | s1$State== 'Igdir'] <- "Northeast Anatolia"
s1$region[s1$State=='Malatya' |
            s1$State== 'Elazig' |
            s1$State== 'Bingol' |
            s1$State== 'Tunceli' |
            s1$State== 'Van' |
            s1$State== 'Mus' |
            s1$State== 'Bitlis' |
            s1$State== 'Hakkari'] <- "Central East Anatolia"
s1$region[s1$State=='Gaziantep' |
            s1$State== 'Adiyaman' |
            s1$State== 'Kilis' |
            s1$State== 'Sanliurfa' |
            s1$State== 'Diyarbakir' |
            s1$State== 'Mardin' |
            s1$State== 'Batman' |
            s1$State== 'Sirnak' |
            s1$State== 'Siirt'] <- "Southeast Anatolia"
table(s1$region, exclude = NULL)

s1$region <- ifelse(is.na(s1$region), "Outside Turkey", s1$region)

#### Supplementary Appendix, Table 16. Study 1 balance testing. ####

library(nnet)

balance_test <-  multinom(factor(Group) ~ beta_sex + beta_age + beta_education + beta_self_left_right + beta_akp + beta_chp + beta_other_right + beta_other_left, data=s1)
summary(balance_test)
lmtest::lrtest(balance_test)
tab_model(balance_test, show.ci = 0.95)

p.adjust(c(0.192, 0.881, 0.972, 0.377, 0.035, 0.041, 0.079, 0.203, 0.619), method = p.adjust.methods, n = length(c(0.192, 0.881, 0.972, 0.377, 0.035, 0.041, 0.079, 0.203, 0.619)))

stargazer(balance_test, type="html", out="s1_balance_test.htm",
          star.cutoffs = c(0.05, 0.01, 0.001),
          digits=2,
          covariate.labels =c("Sex: Male", "Age", "Education", "Left-Right Ideological Placement", "Partisanship: AKP", "Partisanship: CHP", "Partisanship: Other Parties on the Right", "Partisanship: Other Parties on the Left", "Intercept"),
          dep.var.labels   = c("Cancer MS", "Terrorism MS", "Domestic Violence MS"))


#### Supplementary Appendix, Table 19. Study 1 causal mediation analysis models. ####

#### Model 1 ####

out.fit.m01 <- lm(beta_Self_Conservatism ~ beta_Cancer_Prime + beta_Terror_Prime + beta_Domestic_Violence_Prime, data=s1)
summary(out.fit.m01)
tab_model(out.fit.m01, show.ci = 0.95)

confint(out.fit.m01, 'beta_Cancer_Prime', level=0.95)
confint(out.fit.m01, 'beta_Terror_Prime', level=0.95)
confint(out.fit.m01, 'beta_Domestic_Violence_Prime', level=0.95)

C <- diag(4)
model.mc <- glht(out.fit.m01, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

med.fit.m1 <- lm(beta_priming_fear ~ beta_Cancer_Prime + beta_Terror_Prime + beta_Domestic_Violence_Prime, data=s1)
summary(med.fit.m1)
tab_model(med.fit.m1, show.ci = 0.95)

C <- diag(4)
model.mc <- glht(med.fit.m1, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

out.fit.m1 <- lm(beta_Self_Conservatism ~ beta_Cancer_Prime + beta_Terror_Prime + beta_Domestic_Violence_Prime + beta_priming_fear, data=s1)
summary(out.fit.m1)
tab_model(out.fit.m1, show.ci = 0.95)

confint(out.fit.m1, 'beta_Cancer_Prime', level=0.95)
confint(out.fit.m1, 'beta_Terror_Prime', level=0.95)
confint(out.fit.m1, 'beta_Domestic_Violence_Prime', level=0.95)

C <- diag(5)
model.mc <- glht(out.fit.m1, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

library(mediation)

set.seed(123)
med.out.m1.cancer <- mediate(med.fit.m1, out.fit.m1, treat = "beta_Cancer_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.out.m1.cancer)
med.out.m1.cancer$d.avg.ci
 
set.seed(123)
med.out.m1.terror <- mediate(med.fit.m1, out.fit.m1, treat = "beta_Terror_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.out.m1.terror)
med.out.m1.terror$d.avg.ci

set.seed(123)
med.out.m1.domestic.violence <- mediate(med.fit.m1, out.fit.m1, treat = "beta_Domestic_Violence_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.out.m1.domestic.violence)
med.out.m1.domestic.violence$d.avg.ci

#### Model 2 ####

out.fit.m02 <- lm(beta_Self_Conservatism ~ beta_Cancer_Prime + beta_Terror_Prime + beta_Domestic_Violence_Prime + beta_sex + beta_age + beta_education + beta_self_left_right, data=s1)
summary(out.fit.m02)
tab_model(out.fit.m02, show.ci = 0.95)

confint(out.fit.m02, 'beta_Cancer_Prime', level=0.95)
confint(out.fit.m02, 'beta_Terror_Prime', level=0.95)
confint(out.fit.m02, 'beta_Domestic_Violence_Prime', level=0.95)

C <- diag(8)
model.mc <- glht(out.fit.m02, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

med.fit.m2 <- lm(beta_priming_fear ~ beta_Cancer_Prime + beta_Terror_Prime + beta_Domestic_Violence_Prime + beta_sex + beta_age + beta_education + beta_self_left_right, data=s1)
summary(med.fit.m2)
tab_model(med.fit.m2, show.ci = 0.95)

C <- diag(8)
model.mc <- glht(med.fit.m2, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

out.fit.m2 <- lm(beta_Self_Conservatism ~ beta_Cancer_Prime + beta_Terror_Prime + beta_Domestic_Violence_Prime + beta_priming_fear + beta_sex + beta_age + beta_education + beta_self_left_right, data=s1)
summary(out.fit.m2)
tab_model(out.fit.m2, show.ci = 0.95)

confint(out.fit.m2, 'beta_Cancer_Prime', level=0.95)
confint(out.fit.m2, 'beta_Terror_Prime', level=0.95)
confint(out.fit.m2, 'beta_Domestic_Violence_Prime', level=0.95)

C <- diag(9)
model.mc <- glht(out.fit.m2, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

library(mediation)

set.seed(123)
med.out.m2.cancer <- mediate(med.fit.m2, out.fit.m2, treat = "beta_Cancer_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.out.m2.cancer)

set.seed(123)
med.out.m2.terror <- mediate(med.fit.m2, out.fit.m2, treat = "beta_Terror_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.out.m2.terror)

set.seed(123)
med.out.m2.domestic.violence <- mediate(med.fit.m2, out.fit.m2, treat = "beta_Domestic_Violence_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.out.m2.domestic.violence)

library(stargazer)

stargazer(out.fit.m01, med.fit.m1, out.fit.m1, out.fit.m02, med.fit.m2, out.fit.m2, 
          type = "html", title=" ", digits=2, out="s1_model_outputs.htm",
          star.cutoffs = c(0.05, 0.01, 0.001),
          model.numbers = F,
          column.labels = c("Model 0.1.1", "Model 1.1", "Model 1.2", "Model 0.2.1", "Model 2.1", "Model 2.2"),
          covariate.labels = c("Cancer Mortality Salience", 
            "Terrorism Mortality Salience", 
            "Domestic Violence Mortality Salience", 
                               "Mortality Fear",
                               "Sex: Male",
                               "Age", 
                               "Education", 
                               "Left-Right Ideological Placement"),
          add.lines = list(c("Estimation", "OLS", "OLS", "OLS", "OLS", "OLS", "OLS")))

#### Unreported Moderated Mediation by Sex ####
med.fit.m1.sex <- lm(beta_priming_fear ~ beta_sex*beta_Cancer_Prime + beta_sex*beta_Terror_Prime + beta_sex*beta_Domestic_Violence_Prime, data=s1)
summary(med.fit.m1.sex)
tab_model(med.fit.m1.sex, show.ci = 0.95)

C <- diag(8)
model.mc <- glht(med.fit.m1.sex, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

out.fit.m1.sex <- lm(beta_Self_Conservatism ~ beta_sex*beta_Cancer_Prime + beta_sex*beta_Terror_Prime + beta_sex*beta_Domestic_Violence_Prime + beta_sex*beta_priming_fear, data=s1)
summary(out.fit.m1.sex)
tab_model(out.fit.m1.sex, show.ci = 0.95)

C <- diag(10)
model.mc <- glht(out.fit.m1.sex, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

library(mediation)

set.seed(123)
med.fit.m1.cancer.sex <- mediate(med.fit.m1.sex, out.fit.m1.sex, treat = "beta_Cancer_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.fit.m1.cancer.sex)
test.modmed(med.fit.m1.cancer.sex, covariates.1 = list(beta_sex = -0.394485683987275), covariates.2 = list(beta_sex = 0.605514316012725), sims = 10000)

set.seed(123)
med.fit.m1.terror.sex <- mediate(med.fit.m1.sex, out.fit.m1.sex, treat = "beta_Terror_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.fit.m1.terror.sex)
test.modmed(med.fit.m1.terror.sex, covariates.1 = list(beta_sex = -0.394485683987275), covariates.2 = list(beta_sex = 0.605514316012725), sims = 10000)

set.seed(123)
med.fit.m1.domestic.violence.sex <- mediate(med.fit.m1.sex, out.fit.m1.sex, treat = "beta_Domestic_Violence_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.fit.m1.domestic.violence.sex)
test.modmed(med.fit.m1.domestic.violence.sex, covariates.1 = list(beta_sex = -0.394485683987275), covariates.2 = list(beta_sex = 0.605514316012725), sims = 10000)

#### Unreported Moderated Mediation by Age ####
med.fit.m1.age <- lm(beta_priming_fear ~ beta_age*beta_Cancer_Prime + beta_age*beta_Terror_Prime + beta_age*beta_Domestic_Violence_Prime, data=s1)
summary(med.fit.m1.age)
tab_model(med.fit.m1.age, show.ci = 0.95)

C <- diag(8)
model.mc <- glht(med.fit.m1.age, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

out.fit.m1.age <- lm(beta_Self_Conservatism ~ beta_age*beta_Cancer_Prime + beta_age*beta_Terror_Prime + beta_age*beta_Domestic_Violence_Prime + beta_age*beta_priming_fear, data=s1)
summary(out.fit.m1.age)
tab_model(out.fit.m1.age, show.ci = 0.95)

C <- diag(10)
model.mc <- glht(out.fit.m1.age, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

library(mediation)

set.seed(123)
med.fit.m1.cancer.age <- mediate(med.fit.m1.age, out.fit.m1.age, treat = "beta_Cancer_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.fit.m1.cancer.age)
test.modmed(med.fit.m1.cancer.age, covariates.1 = list(beta_age = -0.311738989995775), covariates.2 = list(beta_age = 3.99061440236399), sims = 10000)

set.seed(123)
med.fit.m1.terror.age <- mediate(med.fit.m1.age, out.fit.m1.age, treat = "beta_Terror_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.fit.m1.terror.age)
test.modmed(med.fit.m1.terror.age, covariates.1 = list(beta_age = -0.311738989995775), covariates.2 = list(beta_age = 3.99061440236399), sims = 10000)

set.seed(123)
med.fit.m1.domestic.violence.age <- mediate(med.fit.m1.age, out.fit.m1.age, treat = "beta_Domestic_Violence_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.fit.m1.domestic.violence.age)
test.modmed(med.fit.m1.domestic.violence.age, covariates.1 = list(beta_age = -0.311738989995775), covariates.2 = list(beta_age = 3.99061440236399), sims = 10000)

#### Unreported Moderated Mediation by Left-Right Ideological Placement ####
med.fit.m1.ideology <- lm(beta_priming_fear ~ beta_self_left_right*beta_Cancer_Prime + beta_self_left_right*beta_Terror_Prime + beta_self_left_right*beta_Domestic_Violence_Prime, data=s1)
summary(med.fit.m1.ideology)
tab_model(med.fit.m1.ideology, show.ci = 0.95)

C <- diag(8)
model.mc <- glht(med.fit.m1.ideology, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

out.fit.m1.ideology <- lm(beta_Self_Conservatism ~ beta_self_left_right*beta_Cancer_Prime + beta_self_left_right*beta_Terror_Prime + beta_self_left_right*beta_Domestic_Violence_Prime + beta_self_left_right*beta_priming_fear, data=s1)
summary(out.fit.m1.ideology)
tab_model(out.fit.m1.ideology, show.ci = 0.95)

C <- diag(10)
model.mc <- glht(out.fit.m1.ideology, linfct = C)
summary(model.mc)
BH.OUT <- summary(model.mc, test = adjusted(type = "BH"))
BH.OUT

library(mediation)

set.seed(123)
med.fit.m1.cancer.ideology <- mediate(med.fit.m1.ideology, out.fit.m1.ideology, treat = "beta_Cancer_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.fit.m1.cancer.ideology)
test.modmed(med.fit.m1.cancer.ideology, covariates.1 = list(beta_self_left_right = -0.780801932896792), covariates.2 = list(beta_self_left_right = 0.905351494954022), sims = 10000)

set.seed(123)
med.fit.m1.terror.ideology <- mediate(med.fit.m1.ideology, out.fit.m1.ideology, treat = "beta_Terror_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.fit.m1.terror.ideology)
test.modmed(med.fit.m1.terror.ideology, covariates.1 = list(beta_self_left_right = -0.780801932896792), covariates.2 = list(beta_self_left_right = 0.905351494954022), sims = 10000)

set.seed(123)
med.fit.m1.domestic.violence.ideology <- mediate(med.fit.m1.ideology, out.fit.m1.ideology, treat = "beta_Domestic_Violence_Prime", mediator = "beta_priming_fear", sims = 10000, boot = T, conf.level = 0.95)
summary(med.fit.m1.domestic.violence.ideology)
test.modmed(med.fit.m1.domestic.violence.ideology, covariates.1 = list(beta_self_left_right = -0.780801932896792), covariates.2 = list(beta_self_left_right = 0.905351494954022), sims = 10000)

#### Unreported Serial Mediation Analyses (provided for the interested readers) ####
##### Note to the Readers #####

# I also tested whether the delay tasks (the Big Five Personality Traits) remained stable after the experimental MS and 
# mortality fear and whether alternate causal chains run towards conservatism. 
# Accordingly, I used the serial mediation analysis (MS → mortality fear → delay task (personality traits) → conservative self-identification) 
# by using the R package lavaan (Rosseel et al. 2021). 
# To estimate the power of these SEM models with two serial moderators, I employed a post-hoc Monte Carlo power analysis mc_power_med 
# whose script is written by Schoemann et al. (n.d.).

# Serial mediation analyses, testing the alternative causal chain running from treatments to mortality fear to personality traits to conservatism, 
# revealed that the theorized mechanism was confounded by one unpredicted but the apparently significant and replicable alternate mechanism. 
# The participants under Cancer MS were afraid of cancer mortality, which increased the perception of the self as more neurotic and 
# led to greater conservative self-placement (0.01 standard deviations increase in conservative self-placement with 95% CI [0.003, 0.03] and 82% statistical power).
# It seemed best not to interpret this finding as shifting Neuroticism but rather as considering presenting the self as more or less neurotic 
# in the face of a health threat such as cancer. 
# My cautious interpretation appeared to be supported by the evinced endogeneity between personality traits and politics (see Bakker, Lelkes, and Malka 2021).

# Cited Works

# Bakker, Bert N., Yphtach Lelkes, and Ariel Malka. 2021. ‘Rethinking the link between self-reported personality traits and political preferences’. PsyArXiv. https://doi.org/10.31234/osf.io/a3jf5.
# Rosseel, Yves, Terrence D. Jorgensen, Nicholas Rockwood, Daniel Oberski, Jarrett Byrnes, Leonard Vanbrabant, Victoria Savalei, et al. 2021. Lavaan: Latent variable analysis (version 0.6-9). https://cran.r-project.org/web/packages/lavaan/index.html.
# Schoemann, Alexander M., Aaron J. Boulton, and Stephen D. Short. n.d. Monte Carlo power analysis for indirect effects. https://schoemanna.shinyapps.io/mc_power_med/.

s1$Q20_1[s1$Q20_1==99] <- NA
s1$Q20_2[s1$Q20_2==99] <- NA
s1$Q20_3[s1$Q20_3==99] <- NA
s1$Q20_4[s1$Q20_4==99] <- NA
s1$Q21_1[s1$Q21_1==99] <- NA
s1$Q21_2[s1$Q21_2==99] <- NA
s1$Q21_3[s1$Q21_3==99] <- NA
s1$Q21_4[s1$Q21_4==99] <- NA
s1$Q22_1[s1$Q22_1==99] <- NA
s1$Q22_2[s1$Q22_2==99] <- NA
s1$Q22_3[s1$Q22_3==99] <- NA
s1$Q22_4[s1$Q22_4==99] <- NA
s1$Q23_1[s1$Q23_1==99] <- NA
s1$Q23_2[s1$Q23_2==99] <- NA
s1$Q23_3[s1$Q23_3==99] <- NA
s1$Q23_4[s1$Q23_4==99] <- NA
s1$Q24_1[s1$Q24_1==99] <- NA
s1$Q24_2[s1$Q24_2==99] <- NA
s1$Q24_3[s1$Q24_3==99] <- NA
s1$Q24_4[s1$Q24_4==99] <- NA
s1$Q25_1[s1$Q25_1==99] <- NA
s1$Q25_2[s1$Q25_2==99] <- NA
s1$Q25_3[s1$Q25_3==99] <- NA
s1$Q25_4[s1$Q25_4==99] <- NA
s1$Q26_1[s1$Q26_1==99] <- NA
s1$Q26_2[s1$Q26_2==99] <- NA
s1$Q26_3[s1$Q26_3==99] <- NA
s1$Q26_4[s1$Q26_4==99] <- NA
s1$Q27_1[s1$Q27_1==99] <- NA
s1$Q27_2[s1$Q27_2==99] <- NA
s1$Q27_3[s1$Q27_3==99] <- NA
s1$Q27_4[s1$Q27_4==99] <- NA
s1$Q28_1[s1$Q28_1==99] <- NA
s1$Q28_2[s1$Q28_2==99] <- NA
s1$Q28_3[s1$Q28_3==99] <- NA
s1$Q28_4[s1$Q28_4==99] <- NA
s1$Q29_1[s1$Q29_1==99] <- NA
s1$Q29_2[s1$Q29_2==99] <- NA
s1$Q29_3[s1$Q29_3==99] <- NA
s1$Q29_4[s1$Q29_4==99] <- NA
s1$Q30_1[s1$Q30_1==99] <- NA
s1$Q30_2[s1$Q30_2==99] <- NA
s1$Q30_3[s1$Q30_3==99] <- NA
s1$Q30_4[s1$Q30_4==99] <- NA

names(s1)[names(s1)=="Q20_1"] <- "extraversion_1"
s1$Q21_2R <- 12 - s1$Q21_2
names(s1)[names(s1)=="Q21_2R"] <- "extraversion_2"
names(s1)[names(s1)=="Q22_3"] <- "extraversion_3"
names(s1)[names(s1)=="Q23_4"] <- "extraversion_4"
s1$Q25_1R <- 12 - s1$Q25_1
names(s1)[names(s1)=="Q25_1R"] <- "extraversion_5"
names(s1)[names(s1)=="Q26_2"] <- "extraversion_6"
s1$Q27_3R <- 12 - s1$Q27_3
names(s1)[names(s1)=="Q27_3R"] <- "extraversion_7"
names(s1)[names(s1)=="Q28_4"] <- "extraversion_8"

psych::alpha(s1[s1$Group=="Control", c("extraversion_1", "extraversion_2", "extraversion_3", "extraversion_4",
                   "extraversion_5", "extraversion_6", "extraversion_7", "extraversion_8")])

s1$extraversion <- s1$extraversion_1 + 
  s1$extraversion_2 + 
  s1$extraversion_3 + 
  s1$extraversion_4 + 
  s1$extraversion_5 + 
  s1$extraversion_6 + 
  s1$extraversion_7 + 
  s1$extraversion_8

psych::describeBy(s1$extraversion, s1$Group)

mean(s1$extraversion, na.rm = T)
sd(s1$extraversion, na.rm = T)
s1$beta_extraversion <- (s1$extraversion - mean(s1$extraversion, na.rm = T))/(2*sd(s1$extraversion, na.rm = T))
table(s1$beta_extraversion, exclude = NULL)

s1$Q20_2R <- 12 - s1$Q20_2
names(s1)[names(s1)=="Q20_2R"] <- "agreeableness_1"
names(s1)[names(s1)=="Q21_3"] <- "agreeableness_2"
s1$Q22_4R <- 12 - s1$Q22_4
names(s1)[names(s1)=="Q22_4R"] <- "agreeableness_3"
names(s1)[names(s1)=="Q24_1"] <- "agreeableness_4"
names(s1)[names(s1)=="Q25_2"] <- "agreeableness_5"
s1$Q26_3R <- 12 - s1$Q26_3
names(s1)[names(s1)=="Q26_3R"] <- "agreeableness_6"
names(s1)[names(s1)=="Q27_4"] <- "agreeableness_7"
s1$Q29_1R <- 12 - s1$Q29_1
names(s1)[names(s1)=="Q29_1R"] <- "agreeableness_8"
names(s1)[names(s1)=="Q30_2"] <- "agreeableness_9"

psych::alpha(s1[s1$Group=="Control", c("agreeableness_1", "agreeableness_2", "agreeableness_3", "agreeableness_4",
                   "agreeableness_5", "agreeableness_6", "agreeableness_7", "agreeableness_8", "agreeableness_9")])

s1$agreeableness <- s1$agreeableness_1 + 
  s1$agreeableness_2 + 
  s1$agreeableness_3 + 
  s1$agreeableness_4 + 
  s1$agreeableness_5 + 
  s1$agreeableness_6 + 
  s1$agreeableness_7 + 
  s1$agreeableness_8 +
  s1$agreeableness_9

psych::describeBy(s1$agreeableness, s1$Group)

mean(s1$agreeableness, na.rm = T)
sd(s1$agreeableness, na.rm = T)
s1$beta_agreeableness <- (s1$agreeableness - mean(s1$agreeableness, na.rm = T))/(2*sd(s1$agreeableness, na.rm = T))
table(s1$beta_agreeableness, exclude = NULL)

names(s1)[names(s1)=="Q20_3"] <- "conscientiousness_1"
s1$Q21_4R <- 12 - s1$Q21_4
names(s1)[names(s1)=="Q21_4R"] <- "conscientiousness_2"
names(s1)[names(s1)=="Q23_1"] <- "conscientiousness_3"
s1$Q24_2R <- 12 - s1$Q24_2
names(s1)[names(s1)=="Q24_2R"] <- "conscientiousness_4"
s1$Q25_3R <- 12 - s1$Q25_3
names(s1)[names(s1)=="Q25_3R"] <- "conscientiousness_5"
names(s1)[names(s1)=="Q26_4"] <- "conscientiousness_6"
names(s1)[names(s1)=="Q28_1"] <- "conscientiousness_7"
names(s1)[names(s1)=="Q29_2"] <- "conscientiousness_8"
s1$Q30_3R <- 12 - s1$Q30_3
names(s1)[names(s1)=="Q30_3R"] <- "conscientiousness_9"

psych::alpha(s1[s1$Group=="Control", c("conscientiousness_1", "conscientiousness_2", "conscientiousness_3", "conscientiousness_4",
                   "conscientiousness_5", "conscientiousness_6", "conscientiousness_7", "conscientiousness_8", "conscientiousness_9")])

s1$conscientiousness <- s1$conscientiousness_1 + 
  s1$conscientiousness_2 + 
  s1$conscientiousness_3 + 
  s1$conscientiousness_4 + 
  s1$conscientiousness_5 + 
  s1$conscientiousness_6 + 
  s1$conscientiousness_7 + 
  s1$conscientiousness_8 +
  s1$conscientiousness_9

psych::describeBy(s1$conscientiousness, s1$Group)

mean(s1$conscientiousness, na.rm = T)
sd(s1$conscientiousness, na.rm = T)
s1$beta_conscientiousness <- (s1$conscientiousness - mean(s1$conscientiousness, na.rm = T))/(2*sd(s1$conscientiousness, na.rm = T))
table(s1$beta_conscientiousness, exclude = NULL)

names(s1)[names(s1)=="Q20_4"] <- "neuroticism_1"
s1$Q22_1R <- 12 - s1$Q22_1
names(s1)[names(s1)=="Q22_1R"] <- "neuroticism_2"
names(s1)[names(s1)=="Q23_2"] <- "neuroticism_3"
names(s1)[names(s1)=="Q24_3"] <- "neuroticism_4"
s1$Q25_4R <- 12 - s1$Q25_4
names(s1)[names(s1)=="Q25_4R"] <- "neuroticism_5"
names(s1)[names(s1)=="Q27_1"] <- "neuroticism_6"
s1$Q28_2R <- 12 - s1$Q28_2
names(s1)[names(s1)=="Q28_2R"] <- "neuroticism_7"
names(s1)[names(s1)=="Q29_3"] <- "neuroticism_8"

psych::alpha(s1[s1$Group=="Control", c("neuroticism_1", "neuroticism_2", "neuroticism_3", "neuroticism_4",
                   "neuroticism_5", "neuroticism_6", "neuroticism_7", "neuroticism_8")])

s1$neuroticism <- s1$neuroticism_1 + 
  s1$neuroticism_2 + 
  s1$neuroticism_3 + 
  s1$neuroticism_4 + 
  s1$neuroticism_5 + 
  s1$neuroticism_6 + 
  s1$neuroticism_7 + 
  s1$neuroticism_8

psych::describeBy(s1$neuroticism, s1$Group)

mean(s1$neuroticism, na.rm = T)
sd(s1$neuroticism, na.rm = T)
s1$beta_neuroticism <- (s1$neuroticism - mean(s1$neuroticism, na.rm = T))/(2*sd(s1$neuroticism, na.rm = T))
table(s1$beta_neuroticism, exclude = NULL)

names(s1)[names(s1)=="Q21_1"] <- "openness_1"
names(s1)[names(s1)=="Q22_2"] <- "openness_2"
names(s1)[names(s1)=="Q23_3"] <- "openness_3"
names(s1)[names(s1)=="Q24_4"] <- "openness_4"
names(s1)[names(s1)=="Q26_1"] <- "openness_5"
names(s1)[names(s1)=="Q27_2"] <- "openness_6"
s1$Q28_3R <- 12 - s1$Q28_3
names(s1)[names(s1)=="Q28_3R"] <- "openness_7"
names(s1)[names(s1)=="Q29_4"] <- "openness_8"
s1$Q30_1R <- 12 - s1$Q30_1
names(s1)[names(s1)=="Q30_1R"] <- "openness_9"
names(s1)[names(s1)=="Q30_4"] <- "openness_10"

psych::alpha(s1[s1$Group=="Control", c("openness_1", "openness_2", "openness_3", "openness_4",
                   "openness_5", "openness_6", "openness_7", "openness_8", "openness_9", "openness_10")])

s1$openness <- s1$openness_1 + 
  s1$openness_2 + 
  s1$openness_3 + 
  s1$openness_4 + 
  s1$openness_5 + 
  s1$openness_6 + 
  s1$openness_7 + 
  s1$openness_8 +
  s1$openness_9 +
  s1$openness_10

psych::describeBy(s1$openness, s1$Group)

mean(s1$openness, na.rm = T)
sd(s1$openness, na.rm = T)
s1$beta_openness <- (s1$openness - mean(s1$openness, na.rm = T))/(2*sd(s1$openness, na.rm = T))
table(s1$beta_openness, exclude = NULL)

#### Unreported Study 1 Descriptive Statistics of the Personality Traits. ####

table1(~ openness + conscientiousness + extraversion + agreeableness + neuroticism | Group, data=s1, 
       render.continuous=my.render.cont, render.categorical=my.render.cat, overall=NULL)

library(lavaan)

base.model.1.acme =
  "
  
  beta_priming_fear ~ c_f*beta_Cancer_Prime + t_f*beta_Terror_Prime + d_f*beta_Domestic_Violence_Prime
  
  beta_openness ~ c_b5op*beta_Cancer_Prime + t_b5op*beta_Terror_Prime + d_b5op*beta_Domestic_Violence_Prime + f_b5op*beta_priming_fear + beta_sex + beta_age + beta_education + beta_self_left_right
  beta_conscientiousness ~ c_b5co*beta_Cancer_Prime + t_b5co*beta_Terror_Prime + d_b5co*beta_Domestic_Violence_Prime + f_b5co*beta_priming_fear + beta_sex + beta_age + beta_education + beta_self_left_right
  beta_extraversion ~ c_b5ex*beta_Cancer_Prime + t_b5ex*beta_Terror_Prime + d_b5ex*beta_Domestic_Violence_Prime + f_b5ex*beta_priming_fear + beta_sex + beta_age + beta_education + beta_self_left_right
  beta_agreeableness ~ c_b5ag*beta_Cancer_Prime + t_b5ag*beta_Terror_Prime + d_b5ag*beta_Domestic_Violence_Prime + f_b5ag*beta_priming_fear + beta_sex + beta_age + beta_education + beta_self_left_right
  beta_neuroticism ~ c_b5ne*beta_Cancer_Prime + t_b5ne*beta_Terror_Prime + d_b5ne*beta_Domestic_Violence_Prime + f_b5ne*beta_priming_fear + beta_sex + beta_age + beta_education + beta_self_left_right
  
  beta_Self_Conservatism ~ c_con*beta_Cancer_Prime + t_con*beta_Terror_Prime + d_con*beta_Domestic_Violence_Prime + f_con*beta_priming_fear +
  b5op_con*beta_openness + b5co_con*beta_conscientiousness + b5ex_con*beta_extraversion + b5ag_con*beta_agreeableness + b5ne_con*beta_neuroticism +
  beta_sex + beta_age + beta_education + beta_self_left_right
  
  acme_cancer_fear := c_f*f_con
  acme_terrorism_fear := t_f*f_con
  acme_domestic_violence_fear := d_f*f_con
  
  acme_cancer_fear_b5op := c_f*f_b5op*b5op_con
  acme_cancer_fear_b5co := c_f*f_b5co*b5co_con
  acme_cancer_fear_b5ex := c_f*f_b5ex*b5ex_con
  acme_cancer_fear_b5ag := c_f*f_b5ag*b5ag_con
  acme_cancer_fear_b5ne := c_f*f_b5ne*b5ne_con
  
  acme_terrorism_fear_b5op := t_f*f_b5op*b5op_con
  acme_terrorism_fear_b5co := t_f*f_b5co*b5co_con
  acme_terrorism_fear_b5ex := t_f*f_b5ex*b5ex_con
  acme_terrorism_fear_b5ag := t_f*f_b5ag*b5ag_con
  acme_terrorism_fear_b5ne := t_f*f_b5ne*b5ne_con
  
  acme_domestic_violence_fear_b5op := d_f*f_b5op*b5op_con
  acme_domestic_violence_fear_b5co := d_f*f_b5co*b5co_con
  acme_domestic_violence_fear_b5ex := d_f*f_b5ex*b5ex_con
  acme_domestic_violence_fear_b5ag := d_f*f_b5ag*b5ag_con
  acme_domestic_violence_fear_b5ne := d_f*f_b5ne*b5ne_con
  
  acme_cancer_b5ne := c_b5ne*b5ne_con
  acme_terrorism_b5ne := t_b5ne*b5ne_con
  acme_domestic_violence_b5ne := d_b5ne*b5ne_con
"

fit.base.model.1.acme <-  sem(base.model.1.acme, data=s1, se = "bootstrap", bootstrap = 10000)
summary(fit.base.model.1.acme, fit.measures=TRUE, standardized=TRUE, rsquare=TRUE, nd = 2)
parameterEstimates(fit.base.model.1.acme, rsquare=T, output = "text", level = 0.95)

psych::describe(s1$beta_Cancer_Prime)
psych::describe(s1$beta_Terror_Prime)
psych::describe(s1$beta_Domestic_Violence_Prime)

psych::describe(s1$beta_priming_fear)

psych::describe(s1$beta_agreeableness)
psych::describe(s1$beta_neuroticism)

psych::describe(s1$beta_Self_Conservatism)

#### Supplementary Appendix, Study 1 Sample Representativeness Tests, Tables 1 to 4. ####

table(s1$age, exclude = NULL)
s1$age_cat <- NA
s1$age_cat[s1$age>=18 & s1$age<=24] <- "18-24"
s1$age_cat[s1$age>=25 & s1$age<=34] <- "25-34"
s1$age_cat[s1$age>=35 & s1$age<=49] <- "35-49"
s1$age_cat[s1$age>=50 & s1$age<=64] <- "50-64"
s1$age_cat[s1$age>=65] <- "65 or Over"
table(s1$age_cat, exclude = NULL)
round(prop.table(table(s1$age_cat, exclude = NULL))*100, 1)

s1_age <- c(79.4, 13, 5.8, 1.3, 0.6)
tr_2018_age <- c(15.3, 21.3, 30.1, 21.1, 12.3)

age_s1_tr_2018 <- cbind(s1_age, tr_2018_age)
row.names(age_s1_tr_2018) <- c("18-24","25-34", "35-49", "50-64","65+")

chisq.test(age_s1_tr_2018)
round(chisq.test(age_s1_tr_2018 )$stdres, 2)
round(rcompanion::cramerV(age_s1_tr_2018), 2)

table(s1$sex, exclude = NULL)
s1$sex_cat <- NA
s1$sex_cat[s1$sex==1] <- "Male"
s1$sex_cat[s1$sex==0] <- "Female"
table(s1$sex_cat, exclude = NULL)
round(prop.table(table(s1$sex_cat, exclude = NULL))*100, 1)

s1_gender <- c(39.0, #male
               59.8) #female

tr_2018_gender <- c(49.2,
                    50.8)

gender_s1_tr_2018 <- cbind(s1_gender, tr_2018_gender)
row.names(gender_s1_tr_2018) <- c("Male","Female")

chisq.test(gender_s1_tr_2018)
round(rcompanion::cramerV(gender_s1_tr_2018), 2)

table(s1$education, exclude = NULL)
s1$education_cat <- NA
s1$education_cat[s1$education==2 | s1$education==3] <- "Below High School Graduate"
s1$education_cat[s1$education==4 | s1$education==5] <- "High School Graduate"
s1$education_cat[s1$education==6 | s1$education==7 | s1$education==8] <- "College Graduate or Above"
s1$education_cat[is.na(s1$education)] <- "Unknown"
table(s1$education_cat, exclude = NULL)
s1$education_cat <- factor(s1$education_cat, levels = c("Below High School Graduate", "High School Graduate",
                                                        "College Graduate or Above", "Unknown"))
round(prop.table(table(s1$education_cat, exclude = NULL))*100, 1)

s1_education <- c(2.2, 73.7, 22.8, 1.3)
tr_2018_education <- c(47.9, 24.2, 18, 0.7)

education_s1_tr_2018 <- cbind(s1_education, tr_2018_education)
row.names(education_s1_tr_2018) <- c("Below High School Graduate",
                                     "High School Graduate",
                                     "College Graduate or Above",
                                     "Unknown")

chisq.test(education_s1_tr_2018)
round(chisq.test(education_s1_tr_2018 )$stdres, 2)
round(rcompanion::cramerV(education_s1_tr_2018), 2)

round(prop.table(table(s1$State, exclude = NULL))*100, 1)

s1_geo <- c(0.9,
            0.3,
            0.5,
            0.7,
            0.3,
            0.3,
            1.4,
            6.6,
            0.9,
            0.3,
            0.2,
            0.2,
            0.1,
            0.2,
            0.6,
            0.1,
            0.6,
            0.4,
            0.2,
            0.1,
            0,
            0,
            0,
            0,
            0.3,
            0.5,
            0.2,
            0.2,
            0.1,
            0,
            0.5,
            0.2,
            0.7,
            0.1,
            0.2,
            4.2,
            0.3,
            2.7,
            40.7,
            0.5,
            0.8,
            1.4,
            0.4,
            0.3,
            3.8,
            3.7,
            0.2,
            0.1,
            0,
            0,
            0,
            0,
            0.1,
            0,
            0.3,
            0.4,
            2,
            0.2,
            0.1,
            0.5,
            0.2,
            0,
            0.1,
            2.4,
            0.2,
            10.2,
            1.3,
            0.3,
            0.2,
            0.1,
            0.3,
            0.2,
            0.3,
            0,
            0.1,
            0,
            1,
            0.6,
            0.4,
            0.1,
            0)

tr_2018_geo <- c(1.9,
                 1.4,
                 0.9,
                 1.3,
                 0.7,
                 0.5,
                 1.3,
                 5.7,
                 1.7,
                 0.8,
                 0.5,
                 0.4,
                 0.4,
                 0.3,
                 0.4,
                 0.5,
                 1,
                 1.1,
                 0.7,
                 0.3,
                 0.3,
                 0.4,
                 0.4,
                 0.1,
                 1,
                 1,
                 0.6,
                 0.4,
                 0.2,
                 0.2,
                 1.3,
                 0.5,
                 1.1,
                 0.3,
                 0.4,
                 3.8,
                 0.3,
                 2.4,
                 18.7,
                 1.3,
                 1.8,
                 2.3,
                 0.6,
                 0.6,
                 2.7,
                 3,
                 0.3,
                 0.2,
                 0.5,
                 0.1,
                 0.1,
                 0.3,
                 0.9,
                 0.3,
                 1.9,
                 1.8,
                 2.2,
                 0.8,
                 0.5,
                 0.6,
                 0.2,
                 0.7,
                 0.3,
                 2.7,
                 0.3,
                 6.9,
                 1.7,
                 0.8,
                 0.7,
                 0.4,
                 0.3,
                 0.3,
                 0.3,
                 0.2,
                 0.5,
                 0.8,
                 1.7,
                 1.3,
                 0.7,
                 0.6,
                 0.5)

geo_s1_tr_2018 <- cbind(s1_geo, tr_2018_geo)
row.names(geo_s1_tr_2018) <- c(
  "Manisa",
  "Aydin",
  "Afyon",
  "Mugla",
  "Kutahya",
  "Usak",
  "Denizli",
  "Izmir",
  "Kayseri",
  "Sivas",
  "Aksaray",
  "Nevsehir",
  "Nigde",
  "Kirsehir",
  "Kirikkale",
  "Yozgat",
  "Malatya",
  "Van",
  "Elazig",
  "Hakkari",
  "Bingol",
  "Bitlis",
  "Mus",
  "Tunceli",
  "Ordu",
  "Trabzon",
  "Giresun",
  "Rize",
  "Artvin",
  "Gumushane",
  "Sakarya",
  "Duzce",
  "Eskisehir",
  "Bilecik",
  "Bolu",
  "Bursa",
  "Yalova",
  "Kocaeli",
  "Istanbul",
  "Kahramanmaras",
  "Hatay",
  "Mersin",
  "Isparta",
  "Osmaniye",
  "Adana",
  "Antalya",
  "Burdur",
  "Igdir",
  "Agri",
  "Ardahan",
  "Bayburt",
  "Erzincan",
  "Erzurum",
  "Kars",
  "Sanliurfa",
  "Diyarbakir",
  "Gaziantep",
  "Mardin",
  "Sirnak",
  "Batman",
  "Kilis",
  "Adiyaman",
  "Siirt",
  "Konya",
  "Karaman",
  "Ankara",
  "Samsun",
  "Tokat",
  "Corum",
  "Amasya",
  "Bartin",
  "Karabuk",
  "Sinop",
  "Cankiri",
  "Kastamonu",
  "Zonguldak",
  "Balikesir",
  "Tekirdag",
  "Canakkale",
  "Edirne",
  "Kirklareli")

chisq.test(geo_s1_tr_2018)
round(chisq.test(geo_s1_tr_2018)$stdres, 2)
round(rcompanion::cramerV(geo_s1_tr_2018), 2)

